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We report the first calculations of nuclear properties near the drip-lines using the spherical 
Hartree-Fock-Bogoliubov mean-field theory with a finite-range force supplemented by continuum 
and particle number projection effects. Calculations were carried out in a basis made of the eigen- 
states of a Woods-Saxon potential computed in a box, thereby garanteeing that continuum effects 
were properly taken into account. Projection of the self-consistent solutions on good particle number 
was carried out after variation, and an approximation of the variation after projection result was 
used. We give the position of the drip-lines and examine neutron densities in neutron-rich nuclei. 
We discuss the sensitivity of nuclear observables upon continuum and particle-number restoration 
effects. 

PACS numbers: PACS numbers: 21.10.Dr, 21.10.Gv, 21.10.Re, 21.30.Fe, 21.60.Ev, 21.60.Jz 



O 
=5 

S3 



(N 
> 
O 
OO 
OO 

o 
d 

o 



One of the main challenges in contemporary nuclear 
structure is to describe nuclei far from the valley of stabil- 
ity. New radioactive ion beam facilities are under intense 
development and when in full exploitation will consider- 
ably increase the wealth of experimental data available 
[lj. From a theoretical point of view, neutron- rich or 
-deficient nuclei present two major challenges. Firstly, 
these nuclei are weakly-bound, and it is expected that 
the influence of the continuum of energy will significantly 
alter many nuclear properties • However, the extent to 
which these properties will change due to the presence 
of the continuum is not so clear. Secondly both the nu- 
clear mean-field theory and the nuclear shell model rely 
on effective interactions whose parameters were mostly 
adjusted on data in stable nuclei. The reliability and 
the robustness of these parametrizations under extreme 
isospin ratios must be questioned, and nuclei far from 
stability thus present us with a unique test-ground [l], Q • 

Much effort was devoted in the recent years to the 
study of weakly-bound nuclei, and the development 
of the corresponding appropriate methods, within the 
framework of the non-relativistic mean-field theory with 
Skyrme interactions 0, 0, M , relativistic mean-field @ 
and the nuclear shell model [7|. In mean- field approaches, 
it was recognized early that the stability of very neutron- 
rich nuclei is caused to a great extent by pairing corre- 
lations, and that a fully self-consistent treatment of the 
latter in the Hartree-Fock-Bogoliubov formalism is re- 
quired § . In the case of the Skyrme mean- field the zero- 
range of the interaction leads to the well-known diver- 
gence problem in the pairing channel, and practical cal- 
culations must resort to phenomenological regularization 
procedures [8[. Moreover, approaches based on the effec- 
tive density-functional theory allow terms to be present 
in the functional, that do not originate from a two-body 
Hamiltonian. The price to pay is the occurrence of di- 
vergences when going beyond mean-field and restoring 
the particle number [9|, [ijj] . These difficulties do not ex- 



ist when finite-range effective interactions are used 
However, as of today, the description of weakly-bound 
nuclei in the framework of finite-range mean-field the- 
ory was not available, as existing codes did not take into 
account continuum effects. 

In this paper, we present the first attempt to extend 
the mean-field approach based on finite-range interac- 
tions of the Gogny type to weakly-bound nuclei, by si- 
multaneously taking into account continuum effects and 
the particle number symmetry. In section I we describe 
our method, which relies on the expansion of the mean- 
field solutions on a realistic basis. In section II is shown 
a selection of results benchmarking our method and illus- 
trating the influence of the continuum in weakly-bound 
nuclei calculations. Projection on good particle number 
is discussed in section III. 



I. FINITE-RANGE MEAN-FIELD THEORIES 
WITH COUPLING TO THE CONTINUUM 

Mean-field calculations with finite-range interactions 
are most conveniently carried out in configuration space. 
In the HFB theory quasi-particle wave- functions are thus 
expanded on a basis of well-known functions. In nu- 
clear structure, the harmonic oscillator basis has always 
played a special role, as it produces localized eigenstates 
that are known analytically. However, a well-known defi- 
ciency of this basis is its improper asymptotic behavior. 
While continuum states are spherical waves and therefore 
decay exponentially, HO functions behave like Gaussian 
functions. Unless the number of basis functions is pro- 
hibitively large, it is therefore impossible to describe all 
features of weakly-bound nuclei with such basis states. 

Apart from the techniques already mentioned [!, [3, [f| 
6] , one simple and elegant way to overcome this problem 
is to use the eigenstates of a "realistic" nuclear potential, 
that is, of a potential that tends to at infinity. Such a 
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potential generates both a sequence of strongly localized 
bound-states and of unlocalized positive-energy states. 
In practical calculations, boundary conditions must be 
set. Outgoing wave boundary conditions with wave num- 
ber k assure that all solutions with positive energy are 
taken into account (within the indispensable discretiza- 
tion of the k values in practical numerical calculations), 
including resonant and non-resonant states. However, 
the wave-functions are not square-integrable, and spe- 
cial techniques must be employed to introduce the com- 
pleteness relation and calculate the matrix elements of 
the interaction on such a basis [ll|, [l2|> [IH, [l4[ • On the 
other hand, vanishing boundary conditions, alternatively 
named box boundary conditions, boil down to selecting 
only those positive-energy states that have a node on the 
boundaries of the domain of integration. Consequently, 
not the full continuum is taken into account, but the 
wave-functions are orthonormal. For nuclear structure 
applications and the description of ground-state proper- 
ties the box technique is sufficient and in fact equivalent 
to the full continuum [l5l ]. 

Our model potential was of the Woods-Saxon type, 
which reads in spherical symmetry: 



V(r) 



1 + cxp[(r - r )/a ] 



(1) 



The analytical solutions to the Schrodinger equation 
with this potential are not known apart from s-waves 
states 16]. However, they can be determined easily by 
well-established numerical integration methods. We en- 
close the potential in a sphere of radius Rbox , discretize 
the space with a mesh size h and integrate the radial 
Schrodinger equation numerically: 
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V(r) 



y ni (r) = (2) 



Our tests showed that a mesh size of h ~ 0.1 fm guaran- 
teed the convergence of the subsequent HFB calculation 
to within a few keV. 

The parameters of the potential do not play a major 
role here, as the eigenfunctions are only used as a basis. 
Therefore, the particular form of the Woods-Saxon po- 
tential can only affect the speed at which convergence is 
reached but not the final HFB result. We noticed that in 
order to obtain a fast convergence across the whole mass 
table, the potential should contain a consequent amount 
of bound-states. We therefore adopted the following 
parameters, that correspond to a realistic parametriza- 
tion of the nuclear mean-field of element Z = 126 [l7l |: 
V = -41.619 MeV, r = 1.000 fm, a = 1.693 fm. 

Of more importance is the cut-off of the basis. Like in 
every numerical implementation of the mean-field theory, 
the basis can not be infinite but must be limited artifi- 
cially. In practice, our basis functions are characterized 
by two quantum numbers, £ and n. Contrary to the har- 
monic oscillator basis, where the main shell number N 
provides a convenient and physical cut-off criterion, in 



the case of the most general potential, the existence of 
such a simple criterion can not be guaranteed. Numer- 
ical tests showed that in heavy nuclei, t max = 15 and 
n-max = 18 allowed very accurate calculations from drip- 
line to drip-line to within 10 keV. For the same required 
accuracy, calculations performed with an energy cut-off 
(all basis states with eigenvalues lower than E cut are re- 
tained) reached convergence for E cut of the order of 400 
MeV. About 92% of such a basis is made of positive- 
energy unlocalized continuum states. 

The Gogny interaction is a finite-range force that con- 
tain a central or Brink-Boeker term, a spin-orbit, density- 
dependent and Coulomb terms. The spin-orbit and 
Coulomb terms are the same as for zero-range Skyrme 
forces 18]. Only a few parametrizations of this interac- 
tion exist. All of our calculations were carried out with 
the D1S force [lj|. We will study in a forthcoming paper 
the behavior of the isovector part of the Gogny interac- 
tion in the region of the neutron-drip-linc. 
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FIG. 1: Neutron (left) and proton (right) single-particle ener- 
gies in 208 Pb calculated in the Hartree-Fock approach with the 
D1S Gogny effective interaction (l m ax = 15 and n max = 18). 
The first two columns in each figure correspond to the har- 
monic oscillator basis with N = 24 shells, with either Gauss- 
Hermite ('Analy.') or numerical ('Box.') integrations. The 
column on the right corresponds to the Woods-Saxon basis. 

Figure Q] shows the single-particle spectrum in 208 Pb 
around the Fermi level as function of the type of basis 
used for the calculation. Benchmark calculations were 
performed in the HO basis using standard Gauss-Hermite 
quadrature formulas for the integrations (left column) 
[201 ]. The middle column shows the results when the 
HO basis wave-functions are obtained by integrating the 
Schrodinger equation for a HO potential in a box and us- 
ing numerical integration techniques to compute the ma- 
trix elements of the Gogny interaction. Results obtained 
in the WS basis are shown in the column on the right. 
As required, all three calculations give nearly identical 
results for bound-states, within the numerical accuracy 
of the integration. 
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Let us notice that, as the convergence in the HO basis 
depends on the number of shells and of the oscillator 
length, the convergence in the WS basis depends on the 
cut-offs Umax and n max (or alternatively E cut ), on the 
parameters of the WS potential itself and on the size of 
the box. A detailed study of the convergence properties 
of the WS basis will be presented elsewhere. 



TABLE I: Table of Fermi level and separation energies in 
drip-line nuclei 
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II. DRIP-LINES AND CONTINUUM EFFECTS 

As an application of the consideration of continuum 
effects in practical calculations we present in Table Q] the 
neutron and proton drip-lines calculated in the spherical 
HFB approach for the Woods-Saxon basis. For the listed 
elements, the two-neutrons S2n and two-protons S^p sep- 
aration energies were calculated, as well as the position of 
the Fermi level. The criterion to define the drip-line was 
either a positive separation energy, or a positive Fermi 
level, whichever criterion was met first. The WS basis 
included all states with I < 15 and n < 18 in a box of 
R box = 20 fm. 

Drip-lines can be significantly pushed away by long- 
range correlations. In particular, rare-earth nuclei are 
found to be strongly deformed either in large-scale Gogny 
HFB calculations using the D1S interaction 2l| and in- 
cluding beyond mean-field correlations [12] as well as in 
Skyrme HFB calculations performed in the transformed 
harmonic oscillator basis with the SLy4 parametrization 
in the mean-field channel and the volume delta pairing 
force [23| . Bearing this fact in mind, we nevertheless ob- 
serve systematic differences on the position of the neutron 
drip-line as compared to Skyrme HFB calculations pub- 
lished in [23| where the drip-line almost systematically 
extends further away by 2-4 neutrons. 
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FIG. 2: (color online) Neutron densities in neutron-rich 
132,150,170 g n isotopes compu ted with the finite range Dlb 
Gogny interaction in the HO and Woods-Saxon bases for a 
box of 20 fm. The curve marked Ref. is the HFB calculation 
in coordinate space for 150 Sn with the SkP force [3|. The un- 
physical asymptotic behavior of densities computed in the HO 
basis is clearly visible. Oscillations in the case of 132 Sn are 
caused by the absence of pairing correlations which forbids 
^-couplings in the expansion. 



Observables like densities are very sensitive to contin- 
uum effects and significantly impact many spectroscopic 
or reaction properties of nuclei. A stringent test of our 
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method is therefore its ability to reproduce the long tails 
of neutron densities in neutron-rich nuclei. There have 
been numerous studies of neutron skins in exotic nuclei in 
the past, and we can therefore benchmark our approach. 

As an example we display in figure [2j the neutron den- 
sities of three very neutron-rich Sn isotopes. For small 
and medium values of the nuclear radius the HO and the 
WS bases give the same density profile. However, for 
large radius the differences are significant. As expected, 
the proper asymptotic of the WS basis wave-functions is 
reflected in the long density tails, ft is important to em- 
phasize that the results are nearly identical to coordinate 
space Skyrme HFB calculations of jl| as can be seen from 
the curve marked Ref., which shows the density in 150 Sn 
obtained with the SkP parametrization of the Skyrme in- 
teraction. This provides an independent validation and 
very robust test of our method to incorporate continuum 
effects through a basis. 

We observe oscillations in the densities for closed-shells 
nuclei like 132 Sn. Indeed, in that case the density is sim- 
ply the sum of N well-localized single-particle states, and 
no contribution from discretized continuum states enters 
the expansion. The threshold where these oscillations ap- 
pear gives in fact a measure of the intrinsic numerical ac- 
curacy with which integrations are performed. This phe- 
nomenon may be distinguished from what was reported 
in [24j where oscillations were caused by large discrep- 
ancies between the upper and lower components in the 
HFB equations. 



III. RESTORATION OF THE PARTICLE 
NUMBER SYMMETRY 

As mentioned in the introduction pairing correlations 
play a relevant role in the description of weakly bound 
nuclei. Up to now most descriptions of the drip-line were 
performed in the HFB approach which has two draw- 
backs. Firstly, in the weak pairing regime it usually col- 
lapses to the unpaired solution and secondly, the HFB 
wave function has the right number of particles only on 
the average. To restore this spontaneously broken sym- 
metry, we project the HFB solution onto good particle 
number using the projector: 



1 



2ji 



J<p(N-N) 



(3) 



according to the procedure described in details in Q. 
As is well-known, the variation after projection (VAP) 
method provides a much better approach than the pro- 
jection after variation (PAV) at the cost of a much higher 
numerical effort, especially with the huge basis used in 
the present calculations. In particular, it prevents pair- 
ing correlations to collapse. Nevertheless one can obtain 
an approximate VAP solution by searching for the min- 
imum in a restricted highly pair-correlated variational 
space |25|. This space is generated by HFB solutions 



\(j)(x)) corresponding to the self-consistent diagonaliza- 
tion of the one-body Hamiltonians hjjpB(x), as a func- 
tion of x. Here }ihfb{x) stands for the HFB Hamil- 
tonian obtained by multiplying, at each iteration, the 
pairing field A of the HFB matrix by the constant x, 
i.e., A — > A(x) = xA, while the particle-hole field h 
remains unchanged. To each wave function \<j>(x)) is as- 
sociated an eigenstate of the particle number operator 
$(a;)) = P N \4>{x)) and an energy 



E N {x) 



(*(x)\H\$(x)) 



(4) 



with H the original two-body Hamiltonian. The varia- 
tional principle guarantees that such a procedure yields a 
minimum as function of x, which is an excellent approx- 
imation to the VAP result [2j|. We call this approach 
RVAP (Restricted VAP). 
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FIG. 3: (color online) Total projected energy E N (x) in 42 Mg 
as function of the effective pairing strength factor x. The 
RVAP calculation was performed in the WS basis. 

To investigate the impact of a better treatment of pair- 
ing correlations on exotic nuclei we applied this restricted 
VAP procedure in the neighborhood of the neutron drip- 
line. Figure [3] shows the example of 42 Mg which is un- 
bound against 2 neutrons emission in standard spherical 
HFB or PAV-HFB calculations. The minimum of the 
curve depicted in Fig. [3] corresponds to a stronger pair- 
ing regime than in standard HFB calculations: the total 
gain in energy is only 460 keV but the pairing energy goes 
from -9.22 MeV to -13.83 MeV from x = 1.0 to x = 1.12. 
At the minimum, the two-neutron separation energy cal- 
culated from the RVAP projected energies of 42 Mg and 
40 Mg becomes negative at S% n — —252 keV (compared 
with Sjn = +40 keV in the un-projected case) and the 
nucleus is bound against 2 neutrons emission. This ef- 
fect is usually typical from nuclei near closed-shells where 
only the VAP method can prevent the HFB theory to col- 
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lapse to the unpaired regime. Near the drip-lines, it also 
results in a visible shift of the drip- lines. It is impor- 
tant to realize that this effect can not appear in a PAV 
calculation. 

The restoration of symmetries broken at the mean-field 
level has therefore a significant and experimentally vis- 
ible effect in exotic neutron-rich nuclei. Variation after 
projection on good particle number consists in exploring 
a different variational space than when the projection is 
carried out after the variation. 

In summary, we present the first results of a general- 
ization of HFB calculations with finite-range interactions 
to exotic weakly-bound nuclei. Mean-field solutions are 
expanded on a basis of the eigenstates of a Woods-Saxon 
potential, which possess the proper asymptotic behavior. 
In stable nuclei, this WS basis is benchmarked against 
standard harmonic oscillator calculations and is found 
to give exactly the same results. In neutron-rich nuclei 



we find that the asymptotic of the neutron densities is 
properly reproduced. The role of the restoration of par- 
ticle number before variation is emphasized. It is found 
that this mechanism can stabilize nuclei at the drip-line 
thereby pushing away the latter. 
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